setwd("/mnt/picea/projects/arabidopsis/jhanson/arabidopsis-nutrition-TOR/T0")
suppressPackageStartupMessages({
library(data.table)
library(DESeq2)
library(gplots)
library(here)
library(hyperSpec)
library(limma)
library(LSD)
library(magrittr)
library(matrixStats)
library(parallel)
library(pander)
library(plotly)
library(RColorBrewer)
library(scatterplot3d)
library(tidyverse)
library(tximport)
library(VennDiagram)
library(vsn)
})
suppressPackageStartupMessages({
source("~/Git/UPSCb/UPSCb-common/src/R/featureSelection.R")
source("~/Git/UPSCb/UPSCb-common/src/R/gopher.R")
source("~/Git/UPSCb/UPSCb-common/src/R/plot.multidensity.R")
source("~/Git/UPSCb/UPSCb-common/src/R/volcanoPlot.R")
})
lfc <- 1
FDR <- 0.01
pal <- c(brewer.pal(8,"Dark2"),1)
pal2 <- brewer.pal(9,"Paired") #require package RColorBrewer
cols <- rainbow(17)
hpal <- colorRampPalette(c("blue","white","red"))(100)
mar <- par("mar")
samples <- read.csv("~/Git/UPSCb/projects/arabidopsis-nutrition-TOR/doc/samples3.csv")
samples %<>% filter(!grepl("P11554_1",SciLifeID)) %>%
filter(! SciLifeID %in% c("P13406_101",
"P14066_128",
"P14066_133",
"P13406_102",
"P14066_131")) %>%
filter(! Nutrition == "PKS") %>%
mutate(Nutrition,Nutrition=relevel(Nutrition,"NPS")) %>%
mutate(AZD,AZD=relevel(AZD,"DMSO"))
samples <- samples[order(samples$Timepoint, samples$Nutrition, samples$AZD),]
samples %<>% mutate(Timepoint,Timepoint=factor(paste0("T",Timepoint)))
filenames <- list.files("../Salmon",
recursive = TRUE,
pattern = "quant.sf",
full.names = TRUE)
names(filenames) <- sub("_S.*","",sapply(strsplit(filenames, "/"), .subset, 3))
filenames <- filenames[match(samples$SciLifeID,names(filenames))]
filenames <-filenames[!is.na(names(filenames))]
samples <- samples[match(names(filenames),samples$SciLifeID),]
samples$Conditions <- factor(paste(samples$Timepoint,samples$Nutrition,samples$AZD,sep="_"))
samples$Batch <- factor(substr(samples$SciLifeID,1,8))
Read the expression at the transcript level
tx <- suppressMessages(tximport(files = filenames,
type = "salmon",
txOut = TRUE))
summarise to genes
tx2gene <- data.frame(TXID=rownames(tx$counts),
GENEID=sub("\\.[0-9]+","",rownames(tx$counts)))
gx <- summarizeToGene(tx,tx2gene=tx2gene)
## summarizing abundance
## summarizing counts
## summarizing length
kg <- round(gx$counts)
Sanity check
stopifnot(all(colnames(kg) == samples$SciLifeID))
#dir.create(file.path("analysis_Tom","salmon"),showWarnings=FALSE,recursive=TRUE)
#write.table(kg,file="/mnt/picea/projects/arabidopsis/jhanson/arabidopsis-nutrition-TOR/analysis_Tom/nutrition-unormalised-gene-expression_data.csv")
#save(kg, samples, file = "/mnt/picea/projects/arabidopsis/jhanson/arabidopsis-nutrition-TOR/analysis_Tom/counts.rda")
sel <- rowSums(kg) == 0
sprintf("%s%% (%s) of %s genes are not expressed",
round(sum(sel) * 100/ nrow(kg),digits=1),
sum(sel),
nrow(kg))
## [1] "17.9% (5793) of 32310 genes are not expressed"
For visualization, the data is submitted to a variance stabilization transformation using DESeq2. The dispersion is estimated independently of the sample type
dds <- DESeqDataSetFromMatrix(
countData = kg,
colData = samples,
design = ~ Conditions)
## converting counts to integer mode
dds <- estimateSizeFactors(dds)
vst <- varianceStabilizingTransformation(dds,blind=FALSE)
vsd <- assay(vst)
vsd <- vsd - min(vsd)
Write out
#write.csv(vst,"/mnt/picea/projects/arabidopsis/jhanson/arabidopsis-nutrition-TOR/analysis_Tom/library-size-normalized_variance-stabilized_data_nutrition.csv")
Principal Component Analysis on the normalized data * Establishment of the PCA
conditions1 <- factor(paste(samples$AZD,samples$Nutrition,sep="_"))
pc <- prcomp(t(vsd))
percent <- round(summary(pc)$importance[2,]*100);percent
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8 PC9 PC10 PC11 PC12 PC13 PC14 PC15 PC16
## 43 24 13 6 4 2 1 1 1 1 0 0 0 0 0 0
## PC17 PC18 PC19 PC20 PC21 PC22 PC23 PC24 PC25 PC26 PC27 PC28 PC29 PC30 PC31 PC32
## 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## PC33 PC34 PC35 PC36 PC37 PC38 PC39 PC40 PC41 PC42 PC43 PC44
## 0 0 0 0 0 0 0 0 0 0 0 0
conds <- droplevels(conditions1)
plot(pc$x[,1],
pc$x[,2],
xlab=paste("Comp. 1 (",percent[1],"%)",sep=""),
ylab=paste("Comp. 2 (",percent[2],"%)",sep=""),
col=pal[as.integer(conds)],
pch=c(15,16,17)[as.factor(samples$Timepoint)],
main="All timepoints")
legend("bottomright",pch=19,
col=pal[1:nlevels(conds)],
legend=levels(conds))
legend("topleft",pch=c(15,16,17),
col="black",
legend=c("T0","T6","T24"))
plot(pc$x[,2],
pc$x[,3],
xlab=paste("Comp. 2 (",percent[2],"%)",sep=""),
ylab=paste("Comp. 3 (",percent[3],"%)",sep=""),
col=pal[as.integer(conds)],
pch=c(15,16,17)[as.factor(samples$Timepoint)],
main="All timepoints")
legend("topleft",pch=19,
col=pal[1:nlevels(conds)],
legend=levels(conds))
legend("bottomleft",pch=c(15,16,17),
col="black",
legend=c("T0","T6","T24"))
plot(pc$x[,1],
pc$x[,3],
xlab=paste("Comp. 1 (",percent[1],"%)",sep=""),
ylab=paste("Comp. 3 (",percent[3],"%)",sep=""),
col=pal[as.integer(conds)],
pch=c(15,16,17)[as.factor(samples$Timepoint)],
main="All timepoints")
legend("bottomleft",pch=19,
col=pal[1:nlevels(conds)],
legend=levels(conds))
legend("bottomright",pch=c(15,16,17),
col="black",
legend=c("T0","T6","T24"))
sel <- samples$Timepoint %in% c("T0","T6","T24")
suppressMessages(dds <- DESeqDataSetFromMatrix(
countData = kg[,sel],
colData = samples[sel,],
design = ~ Conditions))
dds <- DESeq(dds)
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
## -- replacing outliers and refitting for 2 genes
## -- DESeq argument 'minReplicatesForReplace' = 7
## -- original counts are preserved in counts(dds)
## estimating dispersions
## fitting model and testing
vst <- varianceStabilizingTransformation(dds,blind=FALSE)
vsd <- assay(vst)
vsd <- vsd - min(vsd)
The contrast by default is the first one (not Intercept)
resultsNames(dds)
## [1] "Intercept" "Conditions_T24_NP_AZD_vs_T0_T0_0"
## [3] "Conditions_T24_NP_DMSO_vs_T0_T0_0" "Conditions_T24_NPS_AZD_vs_T0_T0_0"
## [5] "Conditions_T24_NPS_DMSO_vs_T0_T0_0" "Conditions_T24_NS_AZD_vs_T0_T0_0"
## [7] "Conditions_T24_NS_DMSO_vs_T0_T0_0" "Conditions_T6_NP_AZD_vs_T0_T0_0"
## [9] "Conditions_T6_NP_DMSO_vs_T0_T0_0" "Conditions_T6_NPS_AZD_vs_T0_T0_0"
## [11] "Conditions_T6_NPS_DMSO_vs_T0_T0_0" "Conditions_T6_NS_AZD_vs_T0_T0_0"
## [13] "Conditions_T6_NS_DMSO_vs_T0_T0_0"
res <- results(dds,name = "Conditions_T6_NPS_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NPS_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NPSDMSO <- rownames(res[cutoffs,])
T6NPSDMSO_Low <- rownames(res[cutoff2,])
T6NPSDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 5181 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 3158 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 2023 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T6_NPS_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NPS_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NPSAZD <- rownames(res[cutoffs,])
T6NPSAZD_Low <- rownames(res[cutoff2,])
T6NPSAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 4557 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2619 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1938 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T6_NS_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NS_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NSDMSO <- rownames(res[cutoffs,])
T6NSDMSO_Low <- rownames(res[cutoff2,])
T6NSDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 3402 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 1916 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1486 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T6_NS_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NS_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NSAZD <- rownames(res[cutoffs,])
T6NSAZD_Low <- rownames(res[cutoff2,])
T6NSAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 4804 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2695 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 2109 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T6_NP_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NP_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NPDMSO <- rownames(res[cutoffs,])
T6NPDMSO_Low <- rownames(res[cutoff2,])
T6NPDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 4907 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 3050 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1857 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T6_NP_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T6_NP_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T6NPAZD <- rownames(res[cutoffs,])
T6NPAZD_Low <- rownames(res[cutoff2,])
T6NPAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 3835 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2311 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1524 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NPS_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NPS_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NPSDMSO <- rownames(res[cutoffs,])
T24NPSDMSO_Low <- rownames(res[cutoff2,])
T24NPSDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 6256 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 3062 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 3194 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NPS_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NPS_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NPSAZD <- rownames(res[cutoffs,])
T24NPSAZD_Low <- rownames(res[cutoff2,])
T24NPSAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 4367 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2093 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 2274 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NS_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NS_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NSDMSO <- rownames(res[cutoffs,])
T24NSDMSO_Low <- rownames(res[cutoff2,])
T24NSDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 1794 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 680 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1114 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NS_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NS_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NSAZD <- rownames(res[cutoffs,])
T24NSAZD_Low <- rownames(res[cutoff2,])
T24NSAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 3437 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 1466 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 1971 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NP_DMSO_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NP_DMSO_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NPDMSO <- rownames(res[cutoffs,])
T24NPDMSO_Low <- rownames(res[cutoff2,])
T24NPDMSO_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 5265 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2629 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 2636 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
res <- results(dds,name = "Conditions_T24_NP_AZD_vs_T0_T0_0")
data <- data.frame(rownames(res), res$log2FoldChange, res$padj)
write.table(data, file = "Conditions_T24_NP_AZD_vs_T0_T0_0.csv", sep = ";", row.names = FALSE, dec=",")
cutoffs <- abs(res$log2FoldChange) >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff1 <- res$log2FoldChange >= lfc & ! is.na(res$padj) & res$padj <= FDR
cutoff2 <- res$log2FoldChange <= -lfc & ! is.na(res$padj) & res$padj <= FDR
T24NPAZD <- rownames(res[cutoffs,])
T24NPAZD_Low <- rownames(res[cutoff2,])
T24NPAZD_High <- rownames(res[cutoff1,])
message(sprintf("There are %s genes that are differentially expressed",sum(cutoffs)))
## There are 4803 genes that are differentially expressed
message(sprintf("There are %s genes that are induced",sum(cutoff1)))
## There are 2525 genes that are induced
message(sprintf("There are %s genes that are repressed",sum(cutoff2)))
## There are 2278 genes that are repressed
DESeq2::plotMA(res)
volcanoPlot(res, lfc=1)
hm <- heatmap.2(t(scale(t(vsd[cutoffs,]))),col=hpal,
Colv=FALSE,dendrogram = "row",trace = "none",labRow = FALSE,
distfun = pearson.dist, labCol = paste(samples$Nutrition,samples$AZD,sep="-")[sel],
hclustfun = function(X){hclust(X,method="ward.D2")},margins = c(6,5))
grid.draw(venn.diagram(list(T6NPSDMSO, T6NSDMSO, T6NPDMSO),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSDMSO_Low, T6NSDMSO_Low, T6NPDMSO_Low),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSDMSO_High, T6NSDMSO_High, T6NPDMSO_High),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSAZD, T6NSAZD, T6NPAZD),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSAZD_Low, T6NSAZD_Low, T6NPAZD_Low),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSAZD_High, T6NSAZD_High, T6NPAZD_High),
filename=NULL,
col=pal[1:3],
category.names=c("Full","No P", "No C")))
grid.draw(venn.diagram(list(T6NPSAZD_Low, T6NPSDMSO_Low, T6NPAZD_Low, T6NPDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Suc x AZD", "-Suc")))
grid.draw(venn.diagram(list(T6NPSAZD_High, T6NPSDMSO_High, T6NPAZD_High, T6NPDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Suc x AZD", "-Suc")))
grid.draw(venn.diagram(list(T6NPSAZD_Low, T6NPSDMSO_Low, T6NSAZD_Low, T6NSDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Phos x AZD", "-Phos")))
grid.draw(venn.diagram(list(T6NPSAZD_High, T6NPSDMSO_High, T6NSAZD_High, T6NSDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Phos x AZD", "-Phos")))
grid.draw(venn.diagram(list(T24NPSAZD_Low, T24NPSDMSO_Low, T24NPAZD_Low, T24NPDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Suc x AZD", "-Suc")))
grid.draw(venn.diagram(list(T24NPSAZD_High, T24NPSDMSO_High, T24NPAZD_High, T24NPDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Suc x AZD", "-Suc")))
grid.draw(venn.diagram(list(T24NPSAZD_Low, T24NPSDMSO_Low, T24NSAZD_Low, T24NSDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Phos x AZD", "-Phos")))
grid.draw(venn.diagram(list(T24NPSAZD_High, T24NPSDMSO_High, T24NSAZD_High, T24NSDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("AZD","Ctrl","-Phos x AZD", "-Phos")))
grid.draw(venn.diagram(list(T6NPSAZD_Low, T6NPSDMSO_Low, T24NPAZD_Low, T24NPDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("T6_AZD","T6_Ctrl","T24 -Suc x AZD", "T24 -Suc")))
grid.draw(venn.diagram(list(T6NPSAZD_High, T6NPSDMSO_High, T24NPAZD_High, T24NPDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("T6_AZD","T6_Ctrl","T24 -Suc x AZD", "T24 -Suc")))
grid.draw(venn.diagram(list(T6NSAZD_Low, T6NSDMSO_Low, T24NPAZD_Low, T24NPDMSO_Low),
filename=NULL,
col=pal[1:4],
category.names=c("T6 -Phos. x AZD","T6 -Phos.","T24 -Suc x AZD", "T24 -Suc")))
grid.draw(venn.diagram(list(T6NSAZD_High, T6NSDMSO_High, T24NPAZD_High, T24NPDMSO_High),
filename=NULL,
col=pal[1:4],
category.names=c("T6 -Phos. x AZD","T6 -Phos.","T24 -Suc x AZD", "T24 -Suc")))
write.table(T6NPSDMSO_Low, file = "T6NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NSDMSO_Low, file = "T6NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NPDMSO_Low, file = "T6NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NPSDMSO_High, file = "T6NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T6NSDMSO_High, file = "T6NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T6NPDMSO_High, file = "T6NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T6NPSAZD_Low, file = "T6NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NSAZD_Low, file = "T6NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NPAZD_Low, file = "T6NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T6NPSAZD_High, file = "T6NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(T6NSAZD_High, file = "T6NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(T6NPAZD_High, file = "T6NPAZD_High.csv", sep = ";", row.names = FALSE)
write.table(T24NPSDMSO_Low, file = "T24NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NSDMSO_Low, file = "T24NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NPDMSO_Low, file = "T24NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NPSDMSO_High, file = "T24NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T24NSDMSO_High, file = "T24NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T24NPDMSO_High, file = "T24NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(T24NPSAZD_Low, file = "T24NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NSAZD_Low, file = "T24NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NPAZD_Low, file = "T24NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(T24NPSAZD_High, file = "T24NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(T24NSAZD_High, file = "T24NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(T24NPAZD_High, file = "T24NPAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPSDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana",alpha=2)
## Loading required package: jsonlite
##
## Attaching package: 'jsonlite'
## The following object is masked from 'package:purrr':
##
## flatten
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPSDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NSDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana",alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$opfam, file = "PFAM_T6NSDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPSDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana",alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NSDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NSDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPSAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana",alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPSAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NSAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NSAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPSAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPSAZD_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NSAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NSAZD_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T6NPAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T6NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T6NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T6NPAZD_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPSDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPSDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NSDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NSDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NSDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPDMSO_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPDMSO_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPDMSO_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPSDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPSDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NSDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NSDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NSDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPDMSO_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPDMSO_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPDMSO_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPSAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPSAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NSAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NSAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NSAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPAZD_High,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPAZD_High.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPAZD_High.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPSAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPSAZD_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NSAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NSAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NSAZD_Low.csv", sep = ";", row.names = FALSE)
GO <- gopher(T24NPAZD_Low,background=rownames(vsd)[rowSums(vsd)>0],url="athaliana", alpha=2)
## No enrichments found in task: mapman
write.table(GO$go, file = "GO_T24NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$kegg, file = "KEGG_T24NPAZD_Low.csv", sep = ";", row.names = FALSE)
write.table(GO$pfam, file = "PFAM_T24NPAZD_Low.csv", sep = ";", row.names = FALSE)
## R version 3.6.1 (2019-07-05)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.3 LTS
##
## Matrix products: default
## BLAS/LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] grid parallel stats4 stats graphics grDevices utils
## [8] datasets methods base
##
## other attached packages:
## [1] jsonlite_1.6 vsn_3.54.0
## [3] VennDiagram_1.6.20 futile.logger_1.4.3
## [5] tximport_1.14.0 forcats_0.4.0
## [7] stringr_1.4.0 dplyr_0.8.3
## [9] purrr_0.3.3 readr_1.3.1
## [11] tidyr_1.0.0 tibble_2.1.3
## [13] tidyverse_1.3.0 scatterplot3d_0.3-41
## [15] RColorBrewer_1.1-2 plotly_4.9.1
## [17] pander_0.6.3 magrittr_1.5
## [19] LSD_4.0-0 limma_3.42.0
## [21] hyperSpec_0.99-20180627 ggplot2_3.2.1
## [23] lattice_0.20-38 here_0.1
## [25] gplots_3.0.1.1 DESeq2_1.26.0
## [27] SummarizedExperiment_1.16.0 DelayedArray_0.12.0
## [29] BiocParallel_1.20.0 matrixStats_0.55.0
## [31] Biobase_2.46.0 GenomicRanges_1.38.0
## [33] GenomeInfoDb_1.22.0 IRanges_2.20.1
## [35] S4Vectors_0.24.0 BiocGenerics_0.32.0
## [37] data.table_1.12.6
##
## loaded via a namespace (and not attached):
## [1] colorspace_1.4-1 rprojroot_1.3-2 htmlTable_1.13.2
## [4] XVector_0.26.0 fs_1.3.1 base64enc_0.1-3
## [7] rstudioapi_0.10 affyio_1.56.0 bit64_0.9-7
## [10] AnnotationDbi_1.48.0 lubridate_1.7.4 xml2_1.2.2
## [13] splines_3.6.1 geneplotter_1.64.0 knitr_1.26
## [16] zeallot_0.1.0 Formula_1.2-3 broom_0.5.2
## [19] annotate_1.64.0 cluster_2.1.0 dbplyr_1.4.2
## [22] BiocManager_1.30.10 compiler_3.6.1 httr_1.4.1
## [25] backports_1.1.5 assertthat_0.2.1 Matrix_1.2-18
## [28] lazyeval_0.2.2 cli_1.1.0 formatR_1.7
## [31] acepack_1.4.1 htmltools_0.4.0 tools_3.6.1
## [34] affy_1.64.0 gtable_0.3.0 glue_1.3.1
## [37] GenomeInfoDbData_1.2.2 Rcpp_1.0.3 cellranger_1.1.0
## [40] vctrs_0.2.0 preprocessCore_1.48.0 gdata_2.18.0
## [43] nlme_3.1-142 xfun_0.11 rvest_0.3.5
## [46] testthat_2.3.0 lifecycle_0.1.0 gtools_3.8.1
## [49] XML_3.98-1.20 zlibbioc_1.32.0 scales_1.1.0
## [52] hms_0.5.2 lambda.r_1.2.4 curl_4.2
## [55] yaml_2.2.0 memoise_1.1.0 gridExtra_2.3
## [58] rpart_4.1-15 latticeExtra_0.6-28 stringi_1.4.3
## [61] RSQLite_2.1.2 highr_0.8 genefilter_1.68.0
## [64] checkmate_1.9.4 caTools_1.17.1.2 rlang_0.4.2
## [67] pkgconfig_2.0.3 bitops_1.0-6 evaluate_0.14
## [70] htmlwidgets_1.5.1 bit_1.1-14 tidyselect_0.2.5
## [73] R6_2.4.1 generics_0.0.2 Hmisc_4.3-0
## [76] DBI_1.0.0 pillar_1.4.2 haven_2.2.0
## [79] foreign_0.8-72 withr_2.1.2 survival_3.1-7
## [82] RCurl_1.95-4.12 nnet_7.3-12 modelr_0.1.5
## [85] crayon_1.3.4 futile.options_1.0.1 KernSmooth_2.23-16
## [88] rmarkdown_1.18 locfit_1.5-9.1 readxl_1.3.1
## [91] blob_1.2.0 reprex_0.3.0 digest_0.6.23
## [94] xtable_1.8-4 munsell_0.5.0 viridisLite_0.3.0